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We develop a computationally tractable method for calculating correlation functions of the finite temperature 
trapped Bose gas that includes the effects of s-wave interactions. Our approach uses a classical field method to 
model the low energy modes and treats the high energy modes using a Hartree-Fock description. We present 
results of first and second order correlation functions, in position and momentum space, for an experimentally 
realistic system in the temperature range of O.GTc to l.OTc. We also characterize the spatial coherence length of 
the system. Our theory should be applicable in the critical region where experiments are now able to measure 
first and second order correlations. 
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I. INTRODUCTION 

Recent experimental developments in ultra-cold gases have 
seen an increase in interest in atomic correlation measure- 
ments analogous to the photon correlations observed in the 
landmark experiments of Hanbury -Brown and Twiss fH. In 
those experiments the second order correlation function (in- 
tensity correlations) of light was measured, revealing that pho- 
tons from a thermal light source are bunched. Correlation 
measurements in atomic systems are of significant interest due 
to interaction effects and that both Bose and Fermi particles 
can be investigated. 

The first experiments with atoms by Yasuda et at. |2j used a 
neutral (bosonic) atomic beam and confirmed atom bunching. 
Using ultra-cold Bose gases local high order correlations have 
been inferred from 3-body decay rates |3 4|, and first order 
coherence has been studied using matter wave interference (5] 
13 and Bragg spectroscopy QO. 

More recently there has been spectacular experimental 
progress in the spatially resolved measurement of second or- 
der correlations in both bosonic and fermionic ultra-cold gases 
im [ini [ni El [H [Ml. two general approaches are used to 
make these measurements. One approach involves directly 
counting atoms |l9l[T0][Il|, the other uses absorption imaging 
to measure density fluctuations lfT2l [131 [l4l . The applications 
of these measurements have included Bose and Fermi gases in 
harmonic traps l9|[l0l[ni[Il], and in optical lattices |[T2][T4ll . 

Theoretical work has included the extension of optical def- 
initions of coherence to atomic samples 1 15 , 16], and examin- 
ing how correlation measurements made after the sample has 
expanded from its confinement potential relate to those of the 
original in situ system |[T7l [TSl . Theoretical predictions for 
coiTelations functions of the trapped Bose gas have included 
studies of the ideal case AT! [18] [T9l and interacting gases 
within mean-field descriptions |[T5] . To date, the most com- 
prehensive calculations of second order correlations in the in- 
teracting system used the Hartree Fock Bogoliubov formalism 
II20I . however in that work only local correlations were inves- 
tigated. A non-local extension of that formalism was used to 
calculate the first order coiTelation function of the quasi-2D 
Bose gas ED 119]. 

Generally speaking, it is a considerable challenge to cal- 



culate these correlation functions in a manner that consis- 
tently includes the effect of interactions and harmonic con- 
finement. Near the critical point, where a mean-field analysis 
is no longer valid, there are currently no reliable calculations 
that have been applied to realistic experimental systems. The 
primary purpose of this paper is to develop a classical field 
formalism for calculating non-local coiTelation functions that 
can describe the critical region. The necessity of such a the- 
ory is motivated by two recent experiments: The measurement 
of the critical exponent for the correlation length by the ETH 
Zurich group [22 j ; and the measurement of the second order 
correlation function near the critical point by the Orsay group 



Our formalism, developed in Sec. |ll| is based upon a 
separation of the modes of the system according to their 
mean occupation: we use the Projected Gross Pitaevskii 
Equation (PGPE) to treat the low energy, highly occupied 
modes, and a mean field Hartree Fock approach for the 
high energy, sparsely occupied modes. The PGPE approach 
non-perturbatively includes interactions between low energy 
modes of the gas and is applicable to the critical region. In- 
deed, its predictions for the shifts in critical temperature |23J 
are in good agreement with experimental measurements 1 24] . 
For the higher energy modes of the system, fluctuation effects 
are less important and a mean-field (Hartree Fock) approach 
is sufficient. 



In Sec. Ill we present results of first and second order corre- 
lation functions for experimentally realistic systems, in posi- 
tion and momentum space coordinates. The momentum space 
correlation function approximately corresponds to that mea- 
sured after expansion (e.g. see llTSl ). The approximate range 
of temperatures we investigate is 0.6Tc to l.OTc, where Tc is 
the critical temperature. We also consider a coherence length 
observable and use this to compare the results of our method 
to those of a (pure) mean-field approach. 
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II. FORMALISM 

A. System and correlation functions 

We take our system to be described by the second quantized 
Hamiltonian 

H = Jdr l^pt(r) (^-|^V^ + Vtrap{r)^ ^(r) 

+ ^Uo^\r)¥{r)^{r)^{r)Y (1) 

where ^(r) is the quantum Bose field operator, and Uq = 
Airfi^a/m is the interaction strength, with a the s-wave scat- 
tering length. The harmonic trapping potential is given by 

Vtrapir) = {ujIx^ + + ujIz^) , (2) 

where {uJx,ojy,uJz} are the trap frequencies. 



this situation the measured correlations are proportional to the 
in situ momentum correlations, with the relationship between 
final (observed) position R and in situ momentum given by 
R/icxp = p/m. For this reason we also develop our for- 
malism to calculate momentum space correlations, defined in 
terms of the momentum field operator <&(p), as 

n(p) = ($t(p)$(p)), (8) 
G\p,p') - ($t(p')$(p)), (9) 
G2(p,p') = ($t(p')$t(p)4(p)4(p')). (10) 

Many aspects of the formalism we present for position and 
momentum space are similar, and in what follows we focus 
primarily on giving a detailed derivation for the position case, 
and will only comment on any important differences that arise 
in the momentum case. 



B. Finite temperature formalism 



1. Position space correlation functions 

The density and unnormalised first and second order corre- 
lation functions (e.g. see 1 15 1) are defined as 



n{r) = (^^(r)^(r)), 
G\r,r') = (^^(r')^(r)), 
G2(r,r') = (#^(r')*^(r)*(r)*(r')), 



(3) 
(4) 
(5) 



where the averages are to be interpreted as thermal averages. 
We note that local first order correlations are equal to the den- 
sity, i.e. 7i(r) — G^{r,r). For the purposes of interpret- 
ing particle correlations it is convenient to introduce the nor- 
malised versions of the correlation functions 



G^{r,r') 
n(r')n(r) 



(6) 
(7) 



The first order correlation function describes phase fluctu- 
ations in the system, and as G^(r,r') is the one-body den- 
sity matrix of the system, off-diagonal long range order in this 
quantity is the defining characteristic of Bose-Einstein con- 
densation ll25l . The second order function is a measure of pair 
correlations in the system, for instance the atoms tendency 
to cluster (bunch) or separate (anti-bunch). Recent experi- 
ments in a superfluid Fermi gas have revealed non-local pair- 
ing through measurements of second order correlations ifTSll . 



2. Momentum space correlation functions 



Often the correlation functions for ultra-cold atom systems 
are measured after time of flight expansion for time ioxp- In 



To describe the trapped Bose gas we split the full field op- 
erator into two parts, a classical field (tpc) representing the 
modes which are highly occupied (referred to as the coher- 
ent region), and a quantum field (ipi) representing the sparsely 
populated modes (the incoherent region). The cutoff between 
the two will be taken to be where the occupation in the mode 
is around five particles. 



(11) 



We describe the details of this decomposition further in the 
following subsections. There has been a significant body of 
work on applications of the classical field methods to zero 
and finite temperature properties of Bose gases ll23l l26l l27l 

El ESI Eol E] |32l [331 Ell ISl [SI El IM^ 

These studies show that this sphtting of the field operator can 
be done in a consistent manner and provides an accurate de- 
scription of experiment (in particular see Ref. |23|). 
Substituting expression ( [TT[ ) into the Eqs. gives 



n(r) 
G\ry) 
G^(r,r') 



n,(r) + n,;(r), (12) 
Gl{v,r') + Gliv,r'), (13) 
Gl{r, r') + G^r, r') + 2Gl{r, r')Gj(r, r') 
+n,(r)nc(r') + n»(r')ne(r), (14) 



where 



n,{r) = (^](r)^,(r)) 
Gi(r,r') = (^](r')V,(r)), 
G^(r,r') = {i;]{r')^Ur)^,{r)4'Ar')) 



(15) 
(16) 
(17) 



with j = {i,c} for the incoherent and coherent regions re- 
spectively. 

Several approximations have been made in deriving Eqs. 
([T2])-([14]). We assume that: 
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(i) the coherent and incoherent regions are uncorrelated, so 
that expectations of mixed terms of coherent and inco- 
herent operators can be factorized, e.g. {'ip\'4'*il}c'i>i) — 
{'iplipi){tp*'tpc)- For this to be a good approximation we 
require that our fields are expanded in a basis that ap- 
proximately diagonalizes the problem at the energy cut 
off. For our purposes the single particle basis of the har- 
monic oscillator potential is satisfactory (also see OTI ). 

(ii) we can neglect the averages of single fields. We note 
that in the coherent region this is justified because we 
do not make the symmetry breaking approximation to 
describe the condensate. 

(iii) we can neglect the averages of anomalous fields, e.g. 
(V'cV'c)- Within the classical field approximation the 
anomalous average of the coherent operator is zero, 
and the anomalous expectation of the incoherent field 
is also zero within the Hartree Fock approximation we 
use here. 



C. Coherent Region 

Atoms in the coherent region are treated with the Projected 
Gross Pitaevskii Equation formalism 1 43 J . In this approach 
the coherent field is expanded in a harmonic oscillator basis 



■0c(r,t) 



E 



(18) 



where n represents the quantum numbers for the harmonic os- 
cillator basis states {ipn{r)}, and c„ are time-dependent com- 
plex amplitudes. We consistently define the coherent region 
(C) by restricting the summation to only include oscillator 
states of (single particle) energy less than a predetermined en- 
ergy cut off (Ccut)- The equation of motion for the coherent 
field is the PGPE 

-y^ + Vtrapir) Vc(r,t) 



ih 



dt 



(19) 



2m 

+7'{C/o|^c(r,i)PVc(r,t)} 

where the projection operator {V) formalises our basis set re- 
striction to the coherent region. The primary approximation 
used to arrive at the PGPE is to neglect interaction with the 
incoherent region (see |33 1), which is expected to be well jus- 
tified in the equilibrium case. 



Equation ( 19 1 can be efficiently and accurately simulated 
in the basis set representation using an appropriately chosen 
Gauss-Hermite quadrature (for more details see [44 1). 

The essence of the PGPE formalism is that the field evolu- 



tion according to Eq. ( 19 1 is ergodic and equilibrium proper- 
ties can be found by time averaging solutions. By ignoring the 
coupling to the incoherent modes, the coherent region forms a 
microcanonical system whose equilibrium states are specified 
by the constants of motion, most importantly energy 



2m 



(20) 



and the number of particles 



(21) 



Our simulation procedure is as follows: Initial states are 
prepared with the desired values of Ec and Nc and evolved 
until the sample has thermalized. This is judged by moni- 
toring when the time averages of parameters such as conden- 
sate fraction and temperature have settled down and fluctuate 
about steady state values |50|. The values of /i, T and nc{r) 
which are extracted from the simulations (using methods dis- 
cussed in II27I ') are used further in calculations for the inco- 
herent region. 

Taking Ng samples of the coherent field at times {tj} after 
it has thermalized, we evoke the ergodic hypothesis to evalu- 
ate the correlation functions as 



1 ^= 



(22) 



Gl{r,r') = ^Y^rAr:t,)Mr',t,), (23) 
G^(r,r') = ^Y^\Mr',t,)f\Mr,Wf. (24) 



Correlation functions in momentum space can be immediately 



evaluated using the procedure in ( 22 1-( 24 1, but using the mo- 
mentum field as given by 



P, t) = J dh e-'P-/ Vc(r, t). (25) 



D. Incoherent Region 

Atoms of the incoherent region, i.e. those in high-energy 
sparsely occupied states, are treated with a semiclassical 
phase space approach based on the Hartree Fock meanfield 
analysis. For consistency we only apply this description to the 
volume of phase space complimentary to the coherent region. 

The Wigner distribution for this region is given by 

where fi is the chemical potential and /3 — l/ksT. The en- 
ergy states are given by the Hartree-Fock expression 



E{r,p)^ ^ + Vtrap{r)+2Uo{n,{r) + n,{r)), (27) 
Zm 

and the incoherent density is evaluated as 

n,(r) = J d''pW,{r,p). (28) 

The prime on the integral indicates a restriction in the region 
of phase space coordinates we integrate over, i.e. 



2m 



+ Vtrap(r) > Ecut- 



(29) 
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This condition defines the incoherent region and ensures we 
do not double count contributions already in the coherent re- 
gion. 

The incoherent region formalism is important in determin- 
ing the equilibrium properties of the complete system. It is 
only after self-consistently determining the incoherent density 



using ( 26 1-( 28 1 that we are able to evaluate the total number 
of atoms in the incoherent region as A^i = J d^rni{r), and 
hence the total number of atoms for the system N — Nc + Ni. 

The Wigner function is related to the first order correlation 
function by 1 15| 



Gl{r,r') = J d3pe-'P('-'-')A^^ 



(30) 



Since the phase space coordinates are a diagonal basis for our 
Hartree-Fock treatment of the incoherent region, we can ob- 
tain the second order correlation function as 



(31) 



The momentum space correlation functions for the incoher- 
ent region can be obtained in a similar manner to the position 
space case using the result 

GKp,p') - y"d^re"- (P-P')/'"W, (v, . (32) 

Equation ( [3T] l can also be applied to the momentum space re- 
sult, and we do not repeat it here. 



m. RESULTS 

In this section we present results of the application of our 
formalism to an ultra-cold Bose cloud. The system we con- 
sider consists of approximately 3x10^ ^^Rb atoms confined 
in an anisotropic harmonic trap of frequencies {ujx, cOyjCUz} — 
27r {1, 1, \/8} X 40 s^^. We explore a temperature range of 
approximately O.GTc - l.OTc to investigate the interplay of 
the thermal and condensate clouds, and to see how this af- 
fects the correlation functions. The analytic estimate of the 
critical temperature including finite size and mean-field shifts 
||451 gives Tc ~ 162 nK for this system, however our results 
suggest the actual Tc may be slightly lower, with the system 
attaining ^1% condensate fraction at about 159 nK. 

For our case of a system with external confinement, the 
two-point correlation functions will depend on all coordinates 
Ii5 1 il . Hence a complete characterization of these correlations 
in the 3D system requires six-dimensions. The results we 
present here are for the correlation functions of the full 3D 
system for the case of both points lying on the a;-axis, e.g. 
G^(r = XX, r' = x'x), where x is the unit vector in x direc- 
tion, and this will now be abreviated to just G^{x, x'). The 
geometry of this is indicated in Fig. [T] 



A. Position space correlations 



Figure|2]shows the one body density matrix (G^{x, x')) and 
the normalized first and second order correlation functions for 
temperatures 93 nK [(a)-(c)] and 159 nK [(d)-(f)]. 

In the one body density matrix [Fig. |2](a) and (d)], two 
distinct features are clearly apparent, (i) A narrow ridge runs 
down the diagonal. The peak values of this ridge gives the 
density of the system (recall n{x) — G^(x, x)). We interpret 
the width of the ridge as the length scale over which phase co- 
herence decays for the thermal component of the system, (ii) 
The broad background feature represents the off diagonal long 
range order in the system, which is the defining characteris- 
tic of condensation in interacting systems |25|. For the lower 
temperature result [Fig. |2](a)] the long range order dominates 
due to the large condensate fraction. In contrast, in the higher 
temperature case [Fig. |2](d)] the thermal component is much 
more significant compared to the condensate, which is smaller 
in peak density and spatial extent. This result is close to the 
critical point and the thermal component has a substantial den- 
sity even at the trap center. 

We note that while the background feature arises entirely 
from the coherent region, the ridge has contributions from 
both coherent and incoherent regions. 

In the normalised first order correlations [Fig. |2] (b) and 
(e)], the background and ridge features are still apparent. 
However, normalization emphasizes the thermal component's 




Figure 1: (color online) Instantaneous density slice of classical field 
in xz-plane showing the system and correlation measurement geom- 
etry. The 2;-axis is indicated (white line) and two points x and x' are 
shown. The field and density correlations between these points define 
GI{x, x') and Gc{x, x') respectively. This result is for a system of 
3 X 10" *^Rb atoms in the trap described in the text at a temperature 
of 159 nK. 



contribution at large distances from the trap center. The ridge 
peak value is now unity, as is clear from Eq. (j6]). 

The broadened feature seen in the ridge at large x {x 
x' « IS/irn) is an artifact that arises from the limitations of 
our semiclassical description for the incoherent region. The 
position where this occurs corresponds to the classical turning 
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(a) (d) 




Figure 2: (color online) Position space correlation functions of a harmonically trapped Bose gas. (a)-(c) at 93nK and (d)-(f) 159nK. Other 
parameters: (a)-(c) Nco,id = 153 x 10^ (condensate number), = 1.8 x lO'^, Ni = 1.2 x 10^, e^ut = Mhu^, Ec = Nc x 21.8huj^; (d)-(f) 
Ncond = 3540, iVc = 1.8 X 10", Ni = 2.9 x 10^ e^ut = 36fta;z, Ec = Nc x 20.5ftwz. 



point for the energy cutoff used to define the coherent region 
[|52l. Beyond this point all the modes of the coherent region 
are evanescent, and in examination of Gl{x, x') this manifests 
as an apparent long range order. The local nature of the semi- 
classical approximation means that Gl{x, x') does not cancel 
this feature, as we would expect in a more complete (wave) 
treatment of the incoherent region. With reference to Fig. [2|a) 
and (d) we note that this artifact occurs in a low density region 
and should have a minor effect when results are averaged over 



the whole system. 

We now consider the second order correlations shown in 
Fig. |2](c) and (f), which show a ridge, but no background fea- 
ture. The height of the ridge varies from a value of slightly 
above one up to two. Comparing the second order results 
with the respective first order results, we see that in locations 
where the density is dominated by the condensate, the value 
of g^{x, x) is suppressed from the maximum value of two. 

These observations are consistent with the well known be- 
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Figure 3: (color online) Momentum space correlation functions of a harmonically trapped Bose gas. (a)-(c) at 93nK and (d)-(f) 159nK. 
Other parameters: (a)-(c) Ncond = 153 x 10^, iVe = 1.8 x 10^, Nt = 1.2 x 10^, ecut = 44fttj^, Ec = N^x 21.8ft(.^^; (d)-(f) Ncond = 3540, 
= 1.8 X 10'', Ni = 2.9 X 10^ tcut = 36^^^^, = A^c x 20. 5^^;^. 

havior of g'^{x, x') for photons. For the case of an ideal laser discussed above, is also seen in these results. 
g'^{x,x') — 1, whereas for a thermal light source photon 

bunching occurs with g^{x, x) = 2. ^^^^""^ of correlation functions in momentum space 

We note that the artifact seen in Fig. [2] (b) and (e), and 

I 



Figure b] shows the one body density matrix (G^{px,p'x)) and the normalized first and second order correlation func- 
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tions in momentum space for temperatures 93 nK [(a)-(c)] 
and 159 nK [(d)-(f)]. The dominant feature in G^{px,Px) 
is a sharp spike in momentum which is a signature of con- 
densation. Indeed, we enhance the appearance of the thermal 
component in Fig. |3](a) and (d) by using a logarithmic scale 
in which the density varies by six orders of magnitude. In 
contrast, the position case [Figs. [2|a) and (d)] has a clearly 
discernible thermal component on a linear scale. 

Comparing the low and high temperature results in Figures 
|3ja) and (d) respectively, we observe the peak density of the 
condensate spike to vary by approximately two orders of mag- 
nitude, whereas for the same temperature change the posi- 
tion space density is observed to only change by a factor of 
around two. This observation emphasizes that condensation 
is in some sense a momentum space phenomenon, but also 
marks a rather important difference between the position and 
momentum space correlation functions. The cross-like struc- 
ture extending along the and p'^ axes in Fig. [3](a) is due to 
the enhancement that the large condensate momentum space 
peak provides [this feature is negligible after normalization - 
seeFig.|3](b)]. 

The normalized versions of the correlation functions more 
clearly reveals the non-condensed modes in the system. 
Due to the massive contrast between condensate and non- 
condensate modes in momentum space, experimental mea- 
surement of these correlations will likely prove challenging. 
Indeed, in Ref. |9| results were restricted to above Tc due 
to saturation issues with the detector when a condensate was 
present. 

C. Coherence length 

The coherence length is a measure of first order correlations 
that specifies the typical length scale over which reproducible 
interference fringes might be expected. Recent experiments 
[|22 1 have made detailed measurements of the coherence be- 
tween atoms at different spatial locations within a trapped 
Bose gas using RF-fields to output couple the atoms. This 
approach has the advantage that it avoids volume averaging 
which tends to smear features of the correlation functions (e.g. 
see 1 15 1). This motivates us to consider an on-axis coherence 
length defined as 

^ jdxdx'\G\x,x')\\x~x'f 

2^ dxdx'\G\x,x')\' ^ ' 

which will be used to compare results at different temper- 
atures. This expression is similar in form to the coher- 
ence length defined by Barnett et al. [19J (also see Ref. 
Il46l ). but without volume averaging over the whole system, 
and should be more appropriate for the aforementioned ex- 
periments. For reference, the uniform Boltzmann gas with 
G^{x,x') ~ exp (-7r|a;-x'p/A^^), gives L^, = Xdh/VSn 

where Xdt = h/^/2'KmkT is the thermal de Broglie wave- 
length. 

Figure |4] shows the behaviour of the on-axis coherence 
length for the same system considered in the previous sec- 
tion. We note that this expression for the coherence length is 
dominated by the condensate mode, such that when the con- 
densate fraction is appreciable using the full G^(a;, x') or the 
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50 100 150 200 



T[nK] 

Figure 4: Coherence length L^c calculated using classical field ap- 
proach (circles), and semiclassical HFBP theory (line with dots). For 
reference \db / \/87r is shown as the dashed line. System parameters 
as in Fig. [2] 



condensate contribution G].^^^{x, x') = '^q{x)'^q{x'), where 
^'o(a;) is the condensate mode, yields almost identical results. 
Our results show the coherence length decreases steadily as 
the critical temperature is approached from below. For com- 
parison we have also calculated using a self-consistent 
meanfield approach. Briefly, this approach describes the con- 
densate (for T < Tc) using the Thomas-Fermi approxima- 
tion and the non-condensate using Hartree-Fock Bogoliubov 
Popov (HFBP) theory in the semiclassical limit (e.g. see Ref. 
BSl ). The first order coherence function is then constructed as 
GliFBpi^^x') = Gl^ndi^^x') + Gly,{x,x'), where the ther- 
mal coherence function Gth(a;, x') is obtained via the Wigner 
function as we have explained earlier for the incoherent region 
(see Sec. |llD| but with Ccut ^ 0, also see Ref. L151 ). 

The results of Fig. |4] indicate that the value of calcu- 
lated using meanfield theory is always greater than the value 
obtained from the PGPE simulations for the temperatures 
considered. This primarily occurs because the equilibrium 
condensate fraction in the meanfield results is greater than 
the PGPE case at the same temperature 1 53 1. Nevertheless, 
while there is reasonable agreement between these theories 
for T < Tc, as T ^ Tc their behaviour appears to be rather 
different, with the classical field prediction for decreas- 
ing much less rapidly on the approach to the critical point. 
We note that for the highest temperature classical field re- 
sult (i.e. T — 159 nK) the condensate fraction is ~ 1% and 
the many other low-lying modes will necessarily play an im- 
portant role in the coherence properties of the system in this 
regime. A complete investigation of this behavior is proba- 
bly best done via the correlation length, ^, defined such that 
G^{x, 0) oc x~^ exp{—x/£,), valid for x > X^b, and will be 
the subject of future work. 
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IV. CONCLUSIONS 

We have presented the theoretical development of an effi- 
cient, computationally tractable method for calculating cor- 
relation functions of the finite temperature trapped Bose gas 
in position and momentum space. Our method is robust and 
we have calculated results in a regime comparable to typical 
experiments. 

Our results show the generic characteristics of these func- 
tions, emphasize the striking differences between their behav- 
ior in position and momentum space, and reveal the interplay 
between condensate and thermal components of the system. 
Most current experiments measuring two point correlations 
use expansion imaging, which yields information on the mo- 
mentum space correlations. We have shown that when a con- 
densate is present, this sort of measurement will require sen- 
sitivity to many orders of magnitude variation in atom density 
to yield accurate information on the condensate and thermal 
components. On the other hand, in situ measurements, e.g. 
by output coupling 1221 . may allow direct access to position 
space correlation measurements and require far less sensitiv- 
ity to density variation. 

Having developed a general formalism in this paper, we are 



currently working on applications to several systems of cur- 
rent interest. The primary strength of our method is that the 
low-energy classical field treatment is valid in regimes with 
strong (classical) fluctuations. We are investigating the appli- 
cation of our theory to such fluctuating regimes, motivated by 
two recent experiments: (i) The change of first order coher- 
ence from exponential to algebraic decay has been used as a 
signature of the Berezinskii-Kosterlitz-Thouless transition in 
trapped quasi-2D Bose gas ||6l|42l. (ii) The measurement of 
the critical exponent in the trapped three-dimensional Bose 
Einstein condensation transition by Donner et al. 1221 , yield- 
ing a value of i/ — 0.67 ± 0.13 (compared to the ideal har- 
monically trapped case of = 0.5). 
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